###________________________________________________________________________________________###
###________________________________________________________________________________________###
###________________________________________________________________________________________###
# Script to estimate effects of self monitoring on policy constraint
# Final script updated on November 19, 2018
# Weber and Klar: Final Code to analyze CCAP AND ANES data #
###________________________________________________________________________________________###
###________________________________________________________________________________________###
###________________________________________________________________________________________###

####-------------------------------------------------------------------------------------####
# Code to analyze an GSS Analysis file
####-------------------------------------------------------------------------------------####
rm(list=ls())
require(foreign)
require(dplyr)
require(car)
require(ggplot2)
load("/Users/chrisweber/desktop/Advances Analysis/GSS/data/gss.coded.Rdata")
## Partisan Alignment ##
names(data)
select(data, year, ideo, pid)%>%na.omit()

#Function to compile multiple correlations into a single data frame
corr.data<-OverTimeCorr("pid", "ideo") 
year<-data.frame(year=c(1974:2014))
corr.data<-merge(year, corr.data, by=c("year"), all.x=T)
names(corr.data)<-c("year", "ideo")

corr.data$ideo.fix<-rep(NA, length(corr.data$ideo))
for(i in 1:length(corr.data$ideo)){
  corr.data$ideo.fix[i]<-ifelse(is.na(corr.data$ideo[i]),mean(c(corr.data$ideo[i-1], corr.data$ideo[i+1])), corr.data$ideo[i])  
}
#corr.data$ideo<-ifelse(is.na(corr.data$ideo), 0, corr.data$ideo)


plot1<-ggplot(data = corr.data,
              aes(x = year, 
                  y = ideo))+
  geom_point(position = position_dodge(width = 0.2)) +
  geom_line(data = corr.data,
            aes(x =year, 
                y = ideo.fix, group=1))+
  theme(text=element_text(size=11), 
        axis.text.y=element_text(angle=45))+
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme_bw()+
  # Set the entire chart region to a light gray color
  theme(panel.border=element_rect(colour="white")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.75)) +
  scale_x_continuous(breaks = pretty(corr.data$year, n = 10)) +
  theme(axis.ticks=element_blank()) +
  # Dispose of the legend
  #  theme(legend.position="none") +
  # Set title and axis labels, and format these and tick marks
  ggtitle("Partisan-Ideological Sorting") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  xlab("Year") +
  theme(axis.text.x=element_text(size=11,colour="#535353",face="bold", angle=45)) +
  theme(axis.text.y=element_text(size=11,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Correlation Coefficient", limits=c(0,1))
  


plot1

pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/gss.figure.pdf", width=8, height=7)
plot1
dev.off()

###________________________________________________________________________________________###
###________________________________________________________________________________________###
# CCAP FIGURE #
### CCAP FIGURE ###
rm(list=ls())
load("/Users/chrisweber/Desktop/Advances Analysis/Data/CCAP.data.Rdata")
##############################   CCAP Data ############################## 
sm.data<-subset(sm.data, !is.na(pid.3))  ### Ignore the non-respondents for pid
### Code Affective Polarization ###

sm.data$marry.out<-NA
sm.data$marry.out[sm.data$pid.3=="Democrat"]<-sm.data$marry.rep[sm.data$pid.3=="Democrat"]
sm.data$marry.out[sm.data$pid.3=="Republican"]<-sm.data$marry.dem[sm.data$pid.3=="Republican"]

sm.data$marry.in<-NA
sm.data$marry.in[sm.data$pid.3=="Democrat"]<-sm.data$marry.dem[sm.data$pid.3=="Democrat"]
sm.data$marry.in[sm.data$pid.3=="Republican"]<-sm.data$marry.rep[sm.data$pid.3=="Republican"]

sm.data$move.out<-NA
sm.data$move.out[sm.data$pid.3=="Democrat"]<-sm.data$move.rep[sm.data$pid.3=="Democrat"]
sm.data$move.out[sm.data$pid.3=="Republican"]<-sm.data$move.dem[sm.data$pid.3=="Republican"]

sm.data$move.in<-NA
sm.data$move.in[sm.data$pid.3=="Democrat"]<-sm.data$move.dem[sm.data$pid.3=="Democrat"]
sm.data$move.in[sm.data$pid.3=="Republican"]<-sm.data$move.rep[sm.data$pid.3=="Republican"]

sm.data$feeling.out<-NA
sm.data$feeling.out[sm.data$pid.3=="Democrat"]<-sm.data$rep.feel[sm.data$pid.3=="Democrat"]
sm.data$feeling.out[sm.data$pid.3=="Republican"]<-sm.data$dem.feel[sm.data$pid.3=="Republican"]

sm.data$feeling.in<-NA
sm.data$feeling.in[sm.data$pid.3=="Democrat"]<-sm.data$dem.feel[sm.data$pid.3=="Democrat"]
sm.data$feeling.in[sm.data$pid.3=="Republican"]<-sm.data$rep.feel[sm.data$pid.3=="Republican"]

table(sm.data$move.in)
table(sm.data$move.out)
table(sm.data$move.in)/sum(table(sm.data$move.in))
table(sm.data$move.out)/sum(table(sm.data$move.out))
plot.data<-data.frame(DV=c(sm.data$move.in, sm.data$move.out))
plot.data$Neighbor<-rep(c("Own Party", "Other Party"), each=length(sm.data$move.in))
plot.data<-na.omit(plot.data)

##############################  Neighbor Image  ############################## 
plot1<-  ggplot(plot.data,
                aes(DV))+
  geom_bar(aes(y = ..prop.., fill = factor(..x..)), stat="count") + 
  scale_fill_manual( values=c("#DCDCDC","#DCDCDC","#DCDCDC","#DCDCDC","#DCDCDC","#DCDCDC","#DCDCDC"))+ 
  facet_wrap(~Neighbor)+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Partisan Neighbors") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  theme(axis.text.x=element_text(angle = 90))+
  scale_x_continuous("Emotional Response", breaks=c(1,2,3,4,5,6, 7), labels=c("Happy","2","3","4","5","6","Unhappy"))

pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/plot.social.pdf", width=8, height=7)
plot1
dev.off()


###________________________________________________________________________________________###
###________________________________________________________________________________________###
# ANES Data #
rm(list=ls())
require(foreign)
require(car)
require(nnet)
require(psych)
library(ggplot2)
library(gridExtra)
library(grid)
load("/Users/chrisweber/desktop/Advances Analysis/Data/ANES.data.Rdata")
source("/Users/chrisweber/desktop/Advances Analysis/Analysis/Sm Functions Final.R") ### Source functions
### Generate the sorting measure ###
final.data$sorting<-recode(abs(((final.data$pid*6)+1)-((final.data$ideology*6)+1)), "0=6; 1=5; 2=4; 3=3; 4=2; 6=1")*abs(4-((final.data$pid*6)+1))*abs(4-((final.data$ideology*6)+1))
####Use Post Test Ideology ####
final.data$ideology<-final.data$ideology.post
final.data$int1<-final.data$self.monitoring*final.data$ideology


#########################################################################################################
#########################################################################################################
#########################################################################################################
#########################################################################################################
#########################################################################################################
##### Self Monitoring Moderates the Impact of Symbolic Identities on Policy Outcomes and Party Attachment ####
{
require(pscl)
#Q1: Does self monitoring moderate the impact of ideology #
summary(polr(as.factor(pid)~
           self.monitoring+
           ideology+ideology+
             catholic+no.religion+other.religion+
             female+education+non.white+
             age,
         data=final.data))

m1m<-multinom(pid.3~
                self.monitoring+
                ideology+int1+
                catholic+no.religion+other.religion+
                female+education+non.white+
                age,
                data=final.data)
summary(m1m)

m2m<-lm(affective.polarization~
          self.monitoring+
          ideology+int1+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age,
        data=final.data)

m3m<-lm(affective.polarization2~
          self.monitoring+
          ideology+int1+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age,
        data=final.data)

m4m<-glm(vote.mccain.post~
            self.monitoring+
            ideology+int1+
            catholic+no.religion+other.religion+
            female+education+non.white+
            age,
          data=final.data,
          family=binomial(link="logit"))


voted<-glm(voted~
            self.monitoring+sorting+self.monitoring:sorting+
            catholic+no.religion+other.religion+
            female+education+non.white+
            age,
          data=final.data,
          family=binomial(link="logit"))

psych::alpha(cbind(final.data$efficacy1, final.data$efficacy2, final.data$efficacy3, final.data$efficacy4)) ### Shitty reliability
psych::alpha(cbind(final.data$p1, final.data$p2, final.data$p3, final.data$p4)) ### Shitty reliability
final.data$participation<-rowSums(cbind(final.data$p1, final.data$p2, final.data$p3, final.data$p4, final.data$voted), na.rm=T) ### Shitty reliability

output.coef1<-summary(m1m)$coefficients[1,]
output.coef2<-summary(m1m)$coefficients[2,]
output.se1<-summary(m1m)$standard.errors[1,]
output.se2<-summary(m1m)$standard.errors[2,]


##### Generate the three ANES Table estimates ####
require(dplyr)
output<-cbind(output.coef1, output.se1, output.coef2, output.se2)%>%
  as.data.frame()
detach("package:dplyr")
rownames(output)<-c("Intercept","Self Monitoring", "Ideology",
"Ideology x Self Monitoring", "Catholic", "No Religion", "Other Religion", "Female", "College",
"Non-White", "Age")
names(output)<-c("Democrat", "SE", "Independent",
                    "SE")
require(xtable)
xtable(output, digits=2,
              caption="Party Identification model is a multinomial logistic regression.
              Entries are point estimates and standard errors. 
              Bold entries indicate estimate is two times the standard error. 
              Republican is the excluded category")
dim(model.matrix(m1m))

output<-cbind(summary(m3m)$coefficients[,1:2])
rownames(output)<-c("Intercept","Self Monitoring", "Ideology",
                    "Ideology x Self Monitoring", "Catholic", "No Religion", "Other Religion", "Female", "College",
                    "Non-White", "Age")
xtable(output, digits=2,
       caption="Affective Polarization Models. First model is party polarization (Feelings towards Republicans-Feeling towards Democrats), the second
              is candidate polarization (Feelings towards Mccain-Feeling towards Obama).
              Entries are point estimates and standard errors from OLS regression. 
              Bold entries indicate estimate is two times the standard error. 
              ")

dim(model.matrix(m2m))
dim(model.matrix(m3m))

output<-cbind(summary(m4m)$coefficients[,1:2])
rownames(output)<-c("Intercept","Self Monitoring", "Ideology",
                    "Ideology x Self Monitoring", "Catholic", "No Religion", "Other Religion", "Female", "College",
                    "Non-White", "Age")

xtable(output, digits=2,
       caption="Vote Choice. Logistic regression predicting pr(McCain Vote).
              Entries are point estimates and standard errors from OLS regression. 
              Bold entries indicate estimate is two times the standard error. 
              ")
dim(model.matrix(m4m))

#### Generate the Figures #####
### Figure 1: Party Identification
#########################################################################################################
############################################## Figure 1 (left panel)   ###################################################
#########################################################################################################
#########################################################################################################
detach("package:psych")
require(scales)
####Generate Plots for Two Way Interactions#####
# Figure 1: Party Identification (marginal effects and predictive estimates)
# Figure 2: Voting (Marginal Effects, Predictions)
# Figure 3: Affective polarization
output<-m1m
sm.level<-as.matrix(seq(0,1, 0.1)) ### Simulate across levels of SM
# Predictions
# This code is going to generate three datasets, one for each PID category
data1<-data.frame(matrix(NA, ncol=5))
names(data1)<-c(names(sm.anes.margins(m1m, 0)[[1]]))
data2<-data1
data3<-data1

for(i in 1:length(sm.level)){
  data1<-rbind(data1,sm.anes.margins(m1m, sm.level[i])[[1]])
  data2<-rbind(data2,sm.anes.margins(m1m, sm.level[i])[[2]])
  data3<-rbind(data3,sm.anes.margins(m1m, sm.level[i])[[3]])
}
data<-rbind(data1[-1,], data2[-1,], data3[-1,])
names(data)<-c("min1", "min2", "max2", "max1", "mean")
data$self.monitoring<-rep(sm.level, times=3)
data$party<-rep(c("Democrat", "Independent", "Republican"), each=11)

plot1<-
  ggplot(data=data)+
  geom_ribbon(aes(x=self.monitoring, ymin=min1, ymax=max1, group=party), fill="grey", alpha=0.5)+
  geom_ribbon(aes(x=self.monitoring, ymin=min2, ymax=max2, group=party), fill="grey")+
  geom_line(aes(x=self.monitoring, y=mean, linetype=party))+
  scale_linetype_manual("Party", values=c("twodash", "dotted", "solid"))+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Average Marginal Effect of Ideology") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=13)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Marginal Effect of Ideology", limits=c(-1,1))+
  scale_x_continuous("Self Monitoring")+
  geom_hline(yintercept=0, colour="black")+
  theme(legend.position = "right", legend.background = element_rect(fill="white"))

pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/g1.pdf", width=8, height=7)
plot1
dev.off()

sm.anes.margins(m1m, 0)
sm.anes.margins(m1m, 1)


### Predictions for Voting #####
# Voting
sm.anes.margins2(m4m, 0)
sm.anes.margins2(m4m, 1)

# Marginal effects, across levels of authoritarianism
sm.level<-seq(0,1, 0.1)
# This code is going to generate three datasets, one for each PID category
data1<-data.frame(matrix(NA, ncol=5))
names(data1)<-c(names(sm.anes.margins2(m4m, 0)[[1]]))

for(i in 1:length(sm.level)){
  data1<-rbind(data1,sm.anes.margins2(m4m, sm.level[i])[[1]])
}
data<-data1[-1,]
names(data)<-c("min1", "min2", "max2", "max1", "mean")
data$self.monitoring<-sm.level

plot2<-
  ggplot(data)+
  geom_ribbon(aes(x=self.monitoring, ymin=min1, ymax=max1), fill="grey", alpha=0.5)+
  geom_ribbon(aes(x=self.monitoring, ymin=min2, ymax=max2), fill="grey")+
  geom_line(aes(x=self.monitoring, y=mean))+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Average Marginal Effect of Ideology") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=13)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Marginal Effect of Ideology", limits=c(0,1))+
  scale_x_continuous("Self Monitoring")+
  geom_hline(yintercept=0, colour="black")+
  theme(legend.position = "right", legend.background = element_rect(fill="white"))

pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/g2.pdf", width=8, height=7)
plot2
dev.off()


sm.anes.margins2(m4m, 0)
sm.anes.margins2(m4m, 1)

require(interflex)

inter.kernel(Y = "vote.mccain.post", D = "ideology", X = "self.monitoring", Z = c("catholic", "no.religion", "other.religion", "female",
                                                                                  "education", "non.white", "age"), data = final.data,
             nboots = 50, bw = 1.3, parallel = FALSE, na.rm=T)

inter.binning(Y = "vote.mccain.post", D = "ideology", X = "self.monitoring", Z = c("catholic", "no.religion", "other.religion", "female",
                                                                                  "education", "non.white", "age"), data = final.data,
             na.rm=T)

### Polarization ###
sm.level<-seq(0,1, 0.1)
data1<-data.frame(matrix(NA, ncol=5))
data2<-data.frame(matrix(NA, ncol=5))
names(data1)<-c(names(sm.anes.margins3(m2m, 0)[[1]]))
names(data2)<-c(names(sm.anes.margins3(m2m, 0)[[1]]))

for(i in 1:length(sm.level)){
  data1<-rbind(data1,sm.anes.margins3(m2m, sm.level[i])[[1]])
  data2<-rbind(data2,sm.anes.margins3(m3m, sm.level[i])[[1]])
}
data<-data1[-1,]
names(data)<-c("min1", "min2", "max2", "max1", "mean")
data$self.monitoring<-sm.level
names(data1)<-c("min1", "min2", "max2", "max1", "mean")
names(data2)<-c("min1", "min2", "max2", "max1", "mean")
data1<-data1[-1,]
data2<-data2[-1,]
data1$self.monitoring<-sm.level
data2$self.monitoring<-sm.level


plot3<-
  ggplot(data1)+
  geom_ribbon(aes(x=self.monitoring, ymin=min1, ymax=max1), fill="grey", alpha=0.5)+
  geom_ribbon(aes(x=self.monitoring, ymin=min2, ymax=max2), fill="grey")+
  geom_line(aes(x=self.monitoring, y=mean))+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Average Marginal Effect of Ideology") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=13)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Marginal Effect of Ideology", limits=c(0,1.3))+
  scale_x_continuous("Self Monitoring")+
  geom_hline(yintercept=0, colour="black")+
  theme(legend.position = "right", legend.background = element_rect(fill="white"))

plot4<-
  ggplot(data2)+
  geom_ribbon(aes(x=self.monitoring, ymin=min1, ymax=max1), fill="grey", alpha=0.5)+
  geom_ribbon(aes(x=self.monitoring, ymin=min2, ymax=max2), fill="grey")+
  geom_line(aes(x=self.monitoring, y=mean))+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Average Marginal Effect of Ideology") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=13)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Marginal Effect of Ideology", limits=c(0,1.3))+
  scale_x_continuous("Self Monitoring")+
  geom_hline(yintercept=0, colour="black")+
  theme(legend.position = "right", legend.background = element_rect(fill="white"))


pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/g3.pdf", width=8, height=7)
plot3
dev.off()

pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/g4.pdf", width=8, height=7)
plot4
dev.off()

}
################################################################################################
################################################################################################


#########################################################################################################

## Samara's Data ###

#########################################################################################################

require(readstata13)
data<-read.dta13("/Users/chrisweber/desktop/Advances Analysis/Data/2013 Dataset.dta")

names(data)
table(data$sm1rev)
table(data$sm2)
table(data$sm3rev)
table(data$sm4)
table(data$sm5)
with(data, psych::alpha(cbind(sm1rev, sm3rev, sm4)))
with(data, psych::alpha(cbind(sm2, sm5)))
with(data, psych::alpha(cbind(sm3rev, sm4)))

data$self.monitoring2<-rowMeans(with(data, cbind(sm3rev, sm4)), na.rm=T)
data$self.monitoring2<-rowMeans(with(data, cbind(sm2, sm5)), na.rm=T)
with(data, psych::tetrachoric(cbind(sm1rev, sm2, sm3rev, sm4, sm5)))
with(data, psych::alpha(cbind(sm1rev, sm2, sm3rev, sm4, sm5)))

data$self.monitoring<-data$tsm/5
data$self.monitoring<-rowMeans(cbind(data$sm1rev, data$sm2, data$sm3rev, data$sm4, data$sm5), na.rm=T)
summary(lm(energypriority~self.monitoring+hetero+homo, data))
summary(polr(as.factor(energypriority)~self.monitoring+hetero+homo+hetero:self.monitoring+homo:self.monitoring, data))

#### Analysis of Variance ###
mean(data$energypriority, na.rm=T)
sd(data$energypriority, na.rm=T)

summary(aov(energypriority~as.factor(groupcond), data))
tapply(data$energypriority, data$groupcond, mean, na.rm=T)
tapply(data$energypriority, data$groupcond, sd, na.rm=T)
a<-(lm(energypriority~hetero+homo, data))
summary(a)

#### Generate Figures with Point Estimates ######

a<-(polr(as.factor(energypriority)~self.monitoring+hetero+homo+hetero:self.monitoring+homo:self.monitoring, data))

### Table ####

output<-cbind(summary(a)$coefficients[,1:2])
rownames(output)<-c("Self Monitoring","Heterogeneous Condition", "Homogenous Condition",
                    "Self Monitoring x Heterogeneous", "Self Monitoring x Homogenous", "Cutpoint 1", "Cutpoint 2", "Cutpoint 3", "Cutpoint 4",
                    "Cutpoint 5")

xtable::xtable(output, digits=2,
               caption="Preference for energy policy. Ordered logistic regression results.
               Entries are point estimates and standard errors in parentheses. 
               Bold entries indicate estimate is two times the standard error. 
               ")
dim(model.matrix(a))

#### Generate unique predictions for each group at min and max SM #####
mean.plot<-rbind(
  treatment.effect(a, 0)[[1]],
  treatment.effect(a, 1)[[1]],
  treatment.effect(a, 0)[[2]],
  treatment.effect(a, 1)[[2]],
  treatment.effect(a, 0)[[3]],
  treatment.effect(a, 1)[[3]]
)
names(mean.plot)<-c("min1", "min2", "max2", "max1", "mean")
mean.plot$self.monitoring<-rep(c("low", "high"), times=3)
mean.plot$group<-rep(c("Homogenous", "Heterogeneous", "Control"), each=2)

plot1<-  ggplot(mean.plot,
                aes(x=self.monitoring, y=mean, fill=group)) +
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Treatment Condition", values=c("darkgray", "lightgray", "gray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Predictions across Experimental Conditions") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 1))+
  scale_x_discrete("Self Monitoring")


pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/sam1.pdf", width=8, height=7)
plot1
dev.off()

#### Generate treatment effect plot #####
sm.level<-seq(0,1, 0.1)
plot.data1<-data.frame(matrix(NA, ncol=5))
names(plot.data1)<-c(names(treatment.effect(a, 0)[[4]]))
plot.data2<-data.frame(matrix(NA, ncol=5))
names(plot.data2)<-c(names(treatment.effect(a, 0)[[5]]))

for(i in 1:length(sm.level)){
  plot.data1<-rbind(plot.data1, treatment.effect(a, sm.level[i])[[4]])
  plot.data2<-rbind(plot.data2, treatment.effect(a, sm.level[i])[[5]])
}
plot.data1<-plot.data1[-1,]
plot.data2<-plot.data2[-1,]

plot.data<-rbind(plot.data1, plot.data2)
names(plot.data)<-c("min1", "min2", "max2", "max1", "mean")
plot.data$self.monitoring<-rep(sm.level, times=2)
plot.data$condition<-rep(c("Homogenous", "Hetergeneous"), each=dim(plot.data1)[1])


plot1<-
  ggplot(data=plot.data)+
  geom_ribbon(aes(x=self.monitoring, ymin=min1, ymax=max1, group=condition), fill="grey", alpha=0.5)+
  geom_ribbon(aes(x=self.monitoring, ymin=min2, ymax=max2, group=condition), fill="grey")+
  scale_colour_manual(name="Party", values=c("black", "darkgrey"))+
  geom_line(aes(x=self.monitoring, y=mean, colour=condition))+
  theme_bw() +
  theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Conditional Treatment Effect") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=13)) +
  theme(axis.text.x=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=12,colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=12,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=12,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Treatment Effect", limits=c(-1,1))+
  scale_x_continuous("Self Monitoring")+
  geom_hline(yintercept=0, colour="black")+
  theme(legend.position = "right", legend.background = element_rect(fill="white"))

pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/sam2.pdf", width=8, height=7)
plot1
dev.off()

treatment.effect(a, 1)
treatment.effect(a, 0)

##############################   CCAP Data ############################## 
### CCAP ANALYSIS ###
rm(list=ls())
load("/Users/chrisweber/Desktop/Advances Analysis/Data/CCAP.data.Rdata")
### Load Packages ###
require(psych)
require(car)
require(pscl)
alpha(cbind(sm.data$sm.1, sm.data$sm.2, sm.data$sm.3)) #r=0.68
sm.data$self.monitoring<-rowMeans(cbind(sm.data$sm.1, sm.data$sm.2, sm.data$sm.3), na.rm=T)
sm.data$self.monitoring<-zero.one(sm.data$self.monitoring)
sm.data$republican<-recode(sm.data$pid, "5:7=1; 1:4=0; else=NA")
sm.data$democrat<-recode(sm.data$pid, "1:3=1; 4:7=0; else=NA")
sm.data$independent<-recode(sm.data$pid, "4=1; 1:3=0; 5:7=0; else=NA")
sm.data$mid.1<-recode(sm.data$move.dem, "4=1; 1:3=0; 5:7=0; else=NA")
sm.data$mid.2<-recode(sm.data$move.rep, "4=1; 1:3=0; 5:7=0; else=NA")
sm.data$mid.3<-recode(sm.data$move.neighbor, "4=1; 1:3=0; 5:7=0; else=NA")
sm.data$interaction1<-sm.data$democrat*sm.data$self.monitoring
sm.data$interaction2<-sm.data$independent*sm.data$self.monitoring
### General Effects of Self Monitoring ###

c<-polr(as.factor(marry.rep)~
          self.monitoring+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
d<-polr(as.factor(marry.dem)~
          self.monitoring+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
dim(model.matrix(c))
dim(model.matrix(d))

a<-polr(as.factor(marry.rep)~
          self.monitoring+democrat+independent+interaction1+interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
b<-polr(as.factor(marry.dem)~
          self.monitoring+democrat+independent+interaction1+interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
output.coef1<-summary(a)$coefficients[,1]
output.coef2<-summary(b)$coefficients[,1]
output.se1<-summary(a)$coefficients[,2]
output.se2<-summary(b)$coefficients[,2]

require(dplyr)
output<-cbind(output.coef1, output.se1, output.coef2, output.se2)%>%
  as.data.frame()
detach("package:dplyr")


rownames(output)<-c("Self Monitoring", "Democrat", "Independent",
                    "Self-Monitoring x Democrat", "Self-Monitoring x Independent",
                    "Catholic", "No Religion", "Other Religion", "Female", "College",
                    "Non-White", "Age", "Cutpoint 1", "Cutpoint 2", "Cutpoint 3",
                    "Cutpoint 4", "Cutpoint 5", "Cutpoint 6")
names(output)<-c("Republican Marriage", "SE", "Democratic Marriage",
                 "SE")
require(xtable)
xtable(output, digits=2,
       caption="Party Identification model is a multinomial logistic regression.
       Entries are point estimates and standard errors. 
       Bold entries indicate estimate is two times the standard error. 
       Republican is the excluded category")
### Use this to transition to Partisan Targets
a<-polr(as.factor(marry.rep)~
          self.monitoring+democrat+independent+interaction1+
          interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
b<-polr(as.factor(marry.dem)~
          self.monitoring+democrat+independent+interaction1+
          interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)

dim(model.matrix(a))
dim(model.matrix(b))

#### Plot a Figure Here ####
plot.data<-rbind(
  pred.generator(a,   "ologit", 0, 0)[[2]], ###Happy, DV=Rep, IV=Rep, Low SM
  pred.generator(a,   "ologit", 0, 0)[[1]], ###Happy, DV=Rep, IV=Rep, High SM
  pred.generator(b,   "ologit",  0, 0)[[2]], ###Happy, DV=Dem, IV=Rep, Low SM
  pred.generator(b,   "ologit", 0, 0)[[1]], ###Happy, DV=Dem, IV=Rep, High SM
  
  pred.generator(a,   "ologit", 0, 0)[[4]], ###Unhappy, DV=Rep, IV=Rep, Low SM
  pred.generator(a,   "ologit", 0, 0)[[3]], ###Unhappy, DV=Rep, IV=Rep, High SM
  pred.generator(b,   "ologit", 0, 0)[[4]], ###Unhappy, DV=Dem, IV=Rep, Low SM
  pred.generator(b,   "ologit", 0, 0)[[3]], ###Unhappy, DV=Dem, IV=Rep, High SM
  
  pred.generator(a,   "ologit", 1, 0)[[2]], ###Happy, DV=Rep, IV=Dem, Low SM
  pred.generator(a,   "ologit", 1, 0)[[1]], ###Happy, DV=Rep, IV=Dem, High SM
  pred.generator(b,   "ologit", 1, 0)[[2]], ###Happy, DV=Dem, IV=Dem, Low SM
  pred.generator(b,   "ologit", 1, 0)[[1]], ###Happy, DV=Dem, IV=Dem, High SM
  
  pred.generator(a,   "ologit", 1, 0)[[4]], ###Unhappy, DV=Rep, IV=Dem, Low SM
  pred.generator(a,   "ologit", 1, 0)[[3]], ###Unhappy, DV=Rep, IV=Dem, High SM
  pred.generator(b,   "ologit", 1, 0)[[4]], ###Unhappy, DV=Dem, IV=Dem, Low SM
  pred.generator(b,   "ologit", 1, 0)[[3]], ###Unhappy, DV=Dem, IV=Dem, High SM
  
  pred.generator(a,   "ologit", 0, 1)[[2]], ###Happy, DV=Rep, IV=Ind, Low SM
  pred.generator(a,   "ologit", 0, 1)[[1]], ###Happy, DV=Rep, IV=Ind, High SM
  pred.generator(b,   "ologit", 0, 1)[[2]], ###Happy, DV=Dem, IV=Ind, Low SM
  pred.generator(b,   "ologit", 0, 1)[[1]], ###Happy, DV=Dem, IV=Ind, High SM
  
  pred.generator(a,   "ologit", 0, 1)[[4]], ###Unhappy, DV=Rep, IV=Ind, Low SM
  pred.generator(a,   "ologit", 0, 1)[[3]], ###Unhappy, DV=Rep, IV=Ind, High SM
  pred.generator(b,   "ologit", 0, 1)[[4]], ###Unhappy, DV=Dem, IV=Ind, Low SM
  pred.generator(b,   "ologit", 0, 1)[[3]] ###Unhappy, DV=Dem, IV=Ind, High SM
  
)
names(plot.data)<-c("min1", "min2", "max2", "max1", "mean")
plot.data$self.monitoring<-rep(c("Low SM", "High SM"), times=12)
plot.data$Emotion<-rep(c("Happy", "Unhappy"), each=4, times=3)
plot.data$Neighbor<-rep(c("Republican Marriage", "Democratic Marriage"), each=2, times=6)
plot.data$Self<-rep(c("Republican", "Democrat", "Independent"), each=8)

pred.generator(a,   "ologit", 0, 0)  ## Republican 
pred.generator(b,   "ologit", 0, 0)  ## Democcrat 

pred.generator(a,   "ologit", 1, 0)  ## Rep
pred.generator(b,   "ologit", 1, 0)  ## Dem

#######################################################

### Figure 8 ###

#######################################################

plot1<-  ggplot(subset(plot.data, Self=="Democrat"),
                aes(x=Emotion, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Neighbor)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Marriage. Democratic Respondents") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  scale_x_discrete("Self Monitoring")

#######################################################
plot2<-  ggplot(subset(plot.data, Self=="Republican"),
                aes(x=Emotion, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Neighbor)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Marriage. Republican Respondents") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  scale_x_discrete("Self Monitoring")


plot3<-  ggplot(subset(plot.data, Self=="Independent"),
                aes(x=Emotion, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Neighbor)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Marriage. Independent Respondents") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  scale_x_discrete("Self Monitoring")

pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/figure8a.pdf", width=8, height=7)
plot1
dev.off()
pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/figure8b.pdf", width=8, height=7)
plot2
dev.off()
pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/figure8c.pdf", width=8, height=7)
plot3
dev.off()

# For Text
pred.generator(a,   "ologit", 0, 0) ###Republicans, Republican Marriage
pred.generator(b,   "ologit", 0, 0) ###Republicans, Democratic Marriage
pred.generator(a,   "ologit", 1, 0) ###Democrats, Republican Marriage
pred.generator(b,   "ologit", 1, 0) ###Democrats, Democrat Marriage



#### Neighbor Models ###

a<-polr(as.factor(move.rep)~
          self.monitoring+democrat+independent+interaction1+
          interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
b<-polr(as.factor(move.dem)~
          self.monitoring+democrat+independent+interaction1+
          interaction2+
          catholic+no.religion+other.religion+
          female+education+non.white+
          age, sm.data)
c<-polr(as.factor(context.pid)~
         self.monitoring+democrat+independent+interaction1+
         interaction2+
         catholic+no.religion+other.religion+
         female+education+non.white+
         age, sm.data)


output.coef1<-summary(a)$coefficients[,1]
output.coef2<-summary(b)$coefficients[,1]
output.se1<-summary(a)$coefficients[,2]
output.se2<-summary(b)$coefficients[,2]

require(dplyr)
output<-cbind(output.coef1, output.se1, output.coef2, output.se2)%>%
  as.data.frame()
detach("package:dplyr")


rownames(output)<-c("Self Monitoring", "Democrat", "Independent",
                    "Self-Monitoring x Democrat", "Self-Monitoring x Independent",
                    "Catholic", "No Religion", "Other Religion", "Female", "College",
                    "Non-White", "Age", "Cutpoint 1", "Cutpoint 2", "Cutpoint 3",
                    "Cutpoint 4", "Cutpoint 5", "Cutpoint 6")

names(output)<-c("Republican Neighbor", "SE", "Democratic Neighbor",
                 "SE")


require(xtable)
xtable(output, digits=2,
       caption="Party Identification model is a multinomial logistic regression.
       Entries are point estimates and standard errors. 
       Bold entries indicate estimate is two times the standard error. 
       Republican is the excluded category")


dim(model.matrix(a))
dim(model.matrix(b))

#######################################################

### Figure 2 ### Neighbor ###

#######################################################


plot.data<-rbind(
  pred.generator(a,   "ologit", 0, 0)[[2]], ###Happy, DV=Rep, IV=Rep, Low SM
  pred.generator(a,   "ologit", 0, 0)[[1]], ###Happy, DV=Rep, IV=Rep, High SM
  pred.generator(b,   "ologit",  0, 0)[[2]], ###Happy, DV=Dem, IV=Rep, Low SM
  pred.generator(b,   "ologit", 0, 0)[[1]], ###Happy, DV=Dem, IV=Rep, High SM
  
  pred.generator(a,   "ologit", 0, 0)[[4]], ###Unhappy, DV=Rep, IV=Rep, Low SM
  pred.generator(a,   "ologit", 0, 0)[[3]], ###Unhappy, DV=Rep, IV=Rep, High SM
  pred.generator(b,   "ologit", 0, 0)[[4]], ###Unhappy, DV=Dem, IV=Rep, Low SM
  pred.generator(b,   "ologit", 0, 0)[[3]], ###Unhappy, DV=Dem, IV=Rep, High SM
  
  pred.generator(a,   "ologit", 1, 0)[[2]], ###Happy, DV=Rep, IV=Dem, Low SM
  pred.generator(a,   "ologit", 1, 0)[[1]], ###Happy, DV=Rep, IV=Dem, High SM
  pred.generator(b,   "ologit", 1, 0)[[2]], ###Happy, DV=Dem, IV=Dem, Low SM
  pred.generator(b,   "ologit", 1, 0)[[1]], ###Happy, DV=Dem, IV=Dem, High SM
  
  pred.generator(a,   "ologit", 1, 0)[[4]], ###Unhappy, DV=Rep, IV=Dem, Low SM
  pred.generator(a,   "ologit", 1, 0)[[3]], ###Unhappy, DV=Rep, IV=Dem, High SM
  pred.generator(b,   "ologit", 1, 0)[[4]], ###Unhappy, DV=Dem, IV=Dem, Low SM
  pred.generator(b,   "ologit", 1, 0)[[3]], ###Unhappy, DV=Dem, IV=Dem, High SM
  
  pred.generator(a,   "ologit", 0, 1)[[2]], ###Happy, DV=Rep, IV=Ind, Low SM
  pred.generator(a,   "ologit", 0, 1)[[1]], ###Happy, DV=Rep, IV=Ind, High SM
  pred.generator(b,   "ologit", 0, 1)[[2]], ###Happy, DV=Dem, IV=Ind, Low SM
  pred.generator(b,   "ologit", 0, 1)[[1]], ###Happy, DV=Dem, IV=Ind, High SM
  
  pred.generator(a,   "ologit", 0, 1)[[4]], ###Unhappy, DV=Rep, IV=Ind, Low SM
  pred.generator(a,   "ologit", 0, 1)[[3]], ###Unhappy, DV=Rep, IV=Ind, High SM
  pred.generator(b,   "ologit", 0, 1)[[4]], ###Unhappy, DV=Dem, IV=Ind, Low SM
  pred.generator(b,   "ologit", 0, 1)[[3]] ###Unhappy, DV=Dem, IV=Ind, High SM
  
)
names(plot.data)<-c("min1", "min2", "max2", "max1", "mean")
plot.data$self.monitoring<-rep(c("Low SM", "High SM"), times=12)
plot.data$Emotion<-rep(c("Happy", "Unhappy"), each=4, times=3)
plot.data$Neighbor<-rep(c("Republican Neighbor", "Democratic Neighbor"), each=2, times=6)
plot.data$Self<-rep(c("Republican", "Democrat", "Independent"), each=8)





plot1<-  ggplot(subset(plot.data, Self=="Democrat"),
                aes(x=Emotion, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Neighbor)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Neighbor Democratic Respondents") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  scale_x_discrete("Self Monitoring")

#######################################################
plot2<-  ggplot(subset(plot.data, Self=="Republican"),
                aes(x=Emotion, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Neighbor)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Neighbor Republican Respondents") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 0.8))+
  scale_x_discrete("Self Monitoring")


plot1
plot2

pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/figure8a.pdf", width=8, height=7)
plot1
dev.off()
pdf("/Users/chrisweber/desktop/Advances Analysis/Paper/figure8b.pdf", width=8, height=7)
plot2
dev.off()


# For Text
pred.generator(a,   "ologit", 0, 0) ###Republicans, Republican N
pred.generator(b,   "ologit", 0, 0) ###Republicans, Republican N
pred.generator(a,   "ologit", 1, 0) ###Democrats, Democrat N
pred.generator(b,   "ologit", 1, 0) ###Democrats, Democrat N

##### Reported Friendships #####

plot.data<-rbind(
  pred.generator.friends(c,   "ologit", 0, 0)[[1]], ### DV=Left, IV=Rep, High SM
  pred.generator.friends(c,   "ologit", 0, 0)[[2]], ### DV=Left, IV=Rep, Low SM
  pred.generator.friends(c,   "ologit", 1, 0)[[1]], ### DV=Left, IV=Dem, High SM
  pred.generator.friends(c,   "ologit", 1, 0)[[2]], ### DV=Left, IV=Dem, Low SM
  pred.generator.friends(c,   "ologit", 0, 0)[[3]], ### DV=Right, IV=Rep, High SM
  pred.generator.friends(c,   "ologit", 0, 0)[[4]], ### DV=Right, IV=Rep, Low SM
  pred.generator.friends(c,   "ologit", 1, 0)[[3]], ### DV=Right, IV=Dem, High SM
  pred.generator.friends(c,   "ologit", 1, 0)[[4]]  ### DV=Right, IV=Dem, Low SM
)
names(plot.data)<-c("min1", "min2", "max2", "max1", "mean")
plot.data$self.monitoring<-rep(c("High SM", "Low SM"), times=4)
plot.data$Friends<-rep(c("Democratic Friends", "Republican Friends"),each=4)
plot.data$Self<-rep(c("Republican", "Democrat"), each=2, times=2)

#######################################################
plot1<-  ggplot(plot.data,
                aes(x=Self, y=mean, fill=self.monitoring)) +
  facet_wrap(~ Friends)+
  geom_bar(stat = "identity", position = position_dodge(width = 0.9))+
  geom_errorbar(aes(ymax=max1, ymin=min1), position = position_dodge(width = 0.9), width = 0.05) +
  scale_fill_manual(name="Self Monitoring", values=c("darkgray", "lightgray"))+ 
  theme_bw() +
  #geom_text(aes(label =self.monitoring), 
  #          position = position_dodge(width = 1.5), angle = 0)
  #theme(legend.position="none") +
  theme(panel.background=element_rect(fill="white")) +
  theme(plot.background=element_rect(fill="white")) +
  theme(panel.border=element_rect(colour="#F0F0F0")) +
  # Format the grid
  theme(panel.grid.major=element_line(colour="#D0D0D0",size=.25)) +
  theme(axis.ticks=element_blank())+
  ggtitle("Reported Friends") +
  theme(plot.title=element_text(face="bold",hjust=-.08,vjust=2,colour="#3C3C3C",size=12)) +
  theme(axis.text.x=element_text(size=9,colour="#535353",face="bold")) +
  theme(axis.text.y=element_text(size=9, colour="#535353",face="bold")) +
  theme(axis.title.y=element_text(size=11,colour="#535353",face="bold",vjust=1.5)) +
  theme(axis.title.x=element_text(size=11,colour="#535353",face="bold",vjust=-.5)) +
  scale_y_continuous("Probability", limits=c(0, 1))+
  theme(axis.text.x=element_text(angle = -45, hjust = 0))+
  scale_x_discrete("Respondent PID")

plot1

pred.generator.friends(c,   "ologit", 1, 0)
pred.generator.friends(c,   "ologit", 0, 0)
pred.generator.friends(c,   "ologit", 0, 1)

pdf("/Users/chrisweber/Desktop/Advances Analysis/Paper/friends.pdf", width=8, height=7)
plot1
dev.off()

### Self Monitoring and Power Analysis ####
install.packages("pwr")
require(pwr)
psych::alpha(cbind(scale(sm.data$move.rep), scale(sm.data$marry.rep), scale(abs(5-sm.data$context.pid)))) #0.77
psych::alpha(cbind(scale(sm.data$move.dem), scale(sm.data$marry.dem), scale(abs(sm.data$context.pid)))) #0.72
sm.data$rep.pref<-rowMeans(cbind(scale(sm.data$move.rep), scale(sm.data$marry.rep), scale(abs(5-sm.data$context.pid))), na.rm=T)
sm.data$dem.pref<-rowMeans(cbind(scale(sm.data$move.dem), scale(sm.data$marry.dem), scale(abs(sm.data$context.pid))), na.rm=T)
cor.test(sm.data$rep.pref, sm.data$dem.pref)  #-0.69
sm.data$rep.pref2<-rowMeans(cbind(scale(sm.data$move.rep), scale(sm.data$marry.rep)), na.rm=T)
sm.data$dem.pref2<-rowMeans(cbind(scale(sm.data$move.dem), scale(sm.data$marry.dem)), na.rm=T)
cor.test(sm.data$rep.pref2, sm.data$dem.pref2)  #-0.45

### Power analysis ####
cor.test(sm.data$pid, sm.data$self.monitoring)  #-0.45
t.test(sm.data$self.monitoring[sm.data$pid.3=="Republican"], sm.data$self.monitoring[sm.data$pid.3=="Democrat"])
                                 
#0.23, 0.28, t=-2.89[848.74], p<0.01#
table(sm.data$pid.3[sm.data$self.monitoring>0.5])

p.out <- pwr.p.test(h= ES.h(p1= 0.75, p2 =0.50),
                    sig.level =0.05, 
                    power = 0.80, 
                    alternative= "greater")

cohen.ES(test ="r", size = "medium")
cor.test(sm.data$rep.pref[sm.data$pid.3=="Republican"], sm.data$self.monitoring[sm.data$pid.3=="Republican"])  #-0.07
cor.test(sm.data$dem.pref[sm.data$pid.3=="Republican"], sm.data$self.monitoring[sm.data$pid.3=="Republican"])  #0.03
cor.test(sm.data$rep.pref[sm.data$pid.3=="Democrat"], sm.data$self.monitoring[sm.data$pid.3=="Democrat"])  #0.13
cor.test(sm.data$dem.pref[sm.data$pid.3=="Democrat"], sm.data$self.monitoring[sm.data$pid.3=="Democrat"])  #-0.22
table(sm.data$pid.3[!is.na(sm.data$self.monitoring)]) #377 Republicans
test<-pwr.r.test(r=0.07, sig.level=0.05, power=0.8, alternative="greater")
pwr.r.test(r=0.03, sig.level=0.05, power=0.8, alternative="greater")
pwr.r.test(r=0.07, sig.level=0.05, power=0.8, alternative="greater")

